import numpy
PHICOEF = {
    4: numpy.asarray([
        numpy.asarray([1.5, -2.5, 0., 1.]) * 2/2,
        numpy.asarray([-0.5, 2.5, -4., 2.]) * 2/2,
    ]),
    6: numpy.asarray([
        numpy.asarray([-5., 13., 5., -25., 0., 12.])/12.,
        numpy.asarray([5., -39., 105., -105., 10., 24.])/24.,
        numpy.asarray([-1., 13., -65., 155., -174., 72.])/24.,
    ]),
    8: numpy.asarray([
        numpy.asarray([7., -25., -35., 161., 28., -280., 0., 144.])/144.,
        numpy.asarray([-7., 75., -273., 315., 252., -630., 28., 240.])/240.,
        numpy.asarray([7., -125., 889., -3185., 5908., -4970., 756., 720.])/720.,
        numpy.asarray([-1., 25., -259., 1435., -4564., 8260., -7776., 2880.])/720.,
    ]),
    10: numpy.asarray([
        numpy.asarray([-9., 41., 126., -654., -441., 3129., 324., -5396., 0., 2880.])/2880.,
        numpy.asarray([3., -41., 180., -138., -945., 1911., 690., -3172., 72., 1440.])/1440.,
        numpy.asarray([-9., 205., -1872., 8610., -19845., 15645., 18342., -34540., 3384., 10080.])/10080.,
        numpy.asarray([9., -287., 3870., -28686., 126945., -339423., 523080., -397684., 71856., 40320.])/40320.,
        numpy.asarray([-1., 41., -726., 7266., -45129., 179529., -454544., 700204., -588240., 201600.])/40320.,
    ]),
}
DPHICOEF = {
    4: numpy.asarray([
        numpy.asarray([9, -10, 0])/2.,
        numpy.asarray([-3, 10, -8])/2.,
    ]),
    6: numpy.asarray([
        numpy.asarray([-25, 52, 15, -50, 0])/12.,
        numpy.asarray([25, -156, 315, -210, 10])/24.,
        numpy.asarray([-5, 52, -195, 310, -174])/24.,
    ]),
    8: numpy.asarray([
        numpy.asarray([49, -150, -175, 644, 84, -560, 0])/144.,
        numpy.asarray([-49, 450, -1365, 1260, 756, -1260, 28])/240.,
        numpy.asarray([49, -750, 4445, -12740, 17724, -9940, 756])/720.,
        numpy.asarray([-7, 150, -1295, 5740, -13692, 16520, -7776])/720.,
    ]),
    10: numpy.asarray([
        numpy.asarray([-81, 328, 882, -3924, -2205, 12516, 972, -10792, 0])/2880.,
        numpy.asarray([27, -328, 1260, -828, -4725, 7644, 2070, -6344, 72])/1440.,
        numpy.asarray([-81, 1640, -13104, 51660, -99225, 62580, 55026, -69080, 3384])/10080.,
        numpy.asarray([81, -2296, 27090, -172116, 634725, -1357692, 1569240, -795368, 71856])/40320.,
        numpy.asarray([-9, 328, -5082, 43596, -225645, 718116, -1363632, 1400408, -588240])/40320.,
    ]),
}
divi = {
  4 : [2, 2],
  6 : [12, 24, 24],
  8 : [144, 240, 720, 720],
  10: [2880, 1440, 10080, 40320, 40320]
}
ddivi = {
  4 : [2, 2],
  6 : [12, 24, 24],
  8 : [144, 240, 720, 720],
  10: [2880, 1440, 10080, 40320, 40320]
}

# import sys
# for k, v in PHICOEF.items():
#   print('{')
#   for i in range(v.shape[0]):
#     sys.stdout.write('  {')
#     for j in range(v.shape[1]):
#       sys.stdout.write("%7d. / %5d." % (numpy.round(v[i][j] * divi[k][i]).astype(numpy.int32), divi[k][i]))
#       if j < v.shape[1] - 1:
#         sys.stdout.write(', ')
#     if i < v.shape[0] - 1:
#       print('  },')
#     else:
#       print('  }')
#   print('},')

# for k, v in DPHICOEF.items():
#   print('{')
#   for i in range(v.shape[0]):
#     sys.stdout.write('  {')
#     for j in range(v.shape[1]):
#       sys.stdout.write("%7d. / %5d." % (numpy.round(v[i][j] * ddivi[k][i]).astype(numpy.int32), ddivi[k][i]))
#       if j < v.shape[1] - 1:
#         sys.stdout.write(', ')
#     if i < v.shape[0] - 1:
#       print('  },')
#     else:
#       print('  }')
#   print('},')

half_order_floor = 1
order = 4
phicoef = PHICOEF[order]
dphicoef = DPHICOEF[order]
def compute_phi1(d):
  absd = numpy.abs(d)
  iabsd = numpy.floor(absd).astype(numpy.int64)
  coefs = phicoef[numpy.minimum(iabsd, half_order_floor)]
  print(coefs)
  ret = coefs[:, 0]
  for i in range(1, order):
      ret = ret * absd + coefs[:, i]
  return ret
def compute_dphi1(d):
    absd = numpy.abs(d)
    #iabsd = numpy.asarray(numpy.floor(absd), dtype=numpy.int64)
    iabsd = numpy.floor(absd).astype(numpy.int64)
    coefs = dphicoef[numpy.minimum(iabsd, half_order_floor)]
    print(coefs)
    ret = coefs[:, 0]
    for i in range(1, order - 1):
        ret = ret * absd + coefs[:, i]
    return ret * numpy.sign(d)
print(compute_phi1([0, .5, 1]))
print(compute_dphi1([0, .5, 1]))

print(compute_phi1([1.5, 2, 2.5]))
print(compute_dphi1([1.5, 2, 2.5]))